Here we model an intervention to address malnutrition among ethnic minority groups in Kon Tum, Vietnam. The model follows observations and discussions and some follow up work by partners of the NIFAM (Nutrition Intervention Forecasting and Monitoring) project in 2022 and 2023. A test run of the intervention will be carried out by the Social Policy Ecology research Institute (SPERI) in Kon Tum with smallholder farmers of the Ca Dong ethnic minority group.

Kon Tum garden and cooking class simulation

We use the mcSimulation function to perform a Monte Carlo simulations to estimate model outputs based on provided parameters and a model function. The Monte Carlo simulation generates a set of estimated model outputs based on random input samples, providing a distribution of potential outcomes.

# Source our model
source("KonTum_Garden_Model.R")
## The legacy packages maptools, rgdal, and rgeos, underpinning the sp package,
## which was just loaded, will retire in October 2023.
## Please refer to R-spatial evolution reports for details, especially
## https://r-spatial.org/r/2023/05/15/evolution4.html.
## It may be desirable to make the sf package available;
## package maintainers should consider adding sf to Suggests:.
## The sp package is now running under evolution status 2
##      (status 2 uses the sf package in place of rgdal)
## 
## Attaching package: 'decisionSupport'
## The following objects are masked _by_ '.GlobalEnv':
## 
##     chance_event, discount, vv
# Ensure consistent results with the random number generator
# not for each 'run' of the MC simulation but for 
# consistency each time we call on the simulation 
set.seed(1234) 

garden_simulation_results <- mcSimulation(
  estimate = estimate_read_csv("inputs_kontum_garden.csv"),
  model_function = kontum_garden_function,
  numberOfModelRuns = 1000, #run 1000 times
  functionSyntax = "plainNames"
)

The way we present NPV values can influence decision makers. The same information presented in different ways can even lead to different decisions. Here we derive a plot that compares the decision options as pure NPV values.

source("functions/plot_distributions.R")
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.2     ✔ readr     2.1.4
## ✔ forcats   1.0.0     ✔ stringr   1.5.0
## ✔ ggplot2   3.4.2     ✔ tibble    3.2.1
## ✔ lubridate 1.9.2     ✔ tidyr     1.3.0
## ✔ purrr     1.0.1     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
plot_distributions(mcSimulation_object = garden_simulation_results, 
                                    vars = c("NPV_garden", "NPV_no_garden"),
                                    method = 'hist_simple_overlay', 
                                    base_size = 7, 
                                    x_axis_name = "Comparative NPV outcomes")

Framing the outcomes

Under Prospect Theory the way we present NPV values can influence decision makers - the same information presented in different ways can lead to different decisions. For example, framing a projected NPV gain as a “reduction in potential loss” might make it more attractive to decision makers due to loss aversion.

Here we plot the distribution for the decision and frame the projected NPV gain for the ‘decision’ (distributions for the two options with the NPV values of the no garden option subtracted from those for the garden). If we display this as a “Reduction in potential loss” it might be more attractive to decision makers due to loss aversion, i.e. the party might put more emphasis on avoiding potential losses than on seeking gains. We can frame our results as a strategy that minimizes losses rather than one that maximizes gains.

plot_distributions(mcSimulation_object = garden_simulation_results, 
                                    vars = "decision",
                                    method = 'hist_simple_overlay', 
                                    base_size = 7,  
                                    x_axis_name = "Reduction in potential loss")

Summary of results for the decision

Summary

Use gt_plt_summary() from {gtExtras}

library(gtExtras)
library(svglite)
# Subset the outputs from the mcSimulation function (y) to summarize only on the variables that we want.
mcSimulation_summary <- data.frame(garden_simulation_results$x[2:38], 
                                 garden_simulation_results$y[1:3])

gt_plt_summary(mcSimulation_summary) 
mcSimulation_summary
1000 rows x 40 cols
Column Plot Overview Missing Mean Median SD
discount_rate 3.49.4 0.0% 6.5 6.5 0.9
size_of_garden 28123 0.0% 75.2 75.2 14.9
CV_value 0.010.53 0.0% 0.3 0.3 0.1
inflation_rate 2.912.2 0.0% 7.5 7.5 1.5
if_community_likes 0.020.98 0.0% 0.5 0.5 0.2
if_effective_manage 0.430.82 0.0% 0.6 0.6 0.1
if_garden_yield_enough 0.160.75 0.0% 0.4 0.4 0.1
if_garden_healthy 0.311.00 0.0% 0.7 0.7 0.1
if_effective_training 0.040.98 0.0% 0.5 0.5 0.2
if_offer_green_space 0.290.99 0.0% 0.7 0.7 0.1
if_reduce_polution 0.070.61 0.0% 0.4 0.4 0.1
if_biophysical_good 0.070.66 0.0% 0.4 0.4 0.1
equipment_cost 14130 0.0% 74.9 75.0 15.2
construction_cost 838 0.0% 22.4 22.5 4.6
garden_designing_costs 7.917.4 0.0% 12.5 12.4 1.5
compost_starting 212 0.0% 7.4 7.4 1.5
worm_starting 0.96.2 0.0% 3.5 3.5 0.9
livestock_costs 0.76.9 0.0% 3.5 3.5 0.9
maintaining_labor 1846 0.0% 32.5 32.4 4.7
seed_costs 0.62.4 0.0% 1.5 1.5 0.3
fertilizer 0.52.4 0.0% 1.5 1.5 0.3
plant_protection 0.86.8 0.0% 3.5 3.5 0.9
child_veg_access 313 0.0% 7.5 7.5 1.5
child_healthier_choices 5658 0.0% 290.6 285.5 119.1
women_veg_access 2.612.6 0.0% 7.5 7.5 1.6
women_healthier_choices 4645 0.0% 297.1 298.6 121.3
elderly_veg_access 2.912.6 0.0% 7.5 7.5 1.5
elderly_healthier_choices 10739 0.0% 299.2 292.0 119.3
others_veg_access 2.512.1 0.0% 7.5 7.5 1.6
others_healthier_choices 6771 0.0% 303.3 301.0 119.7
green_space_value 55234 0.0% 149.9 149.1 30.8
reduce_polution_value 726 0.0% 15.2 15.2 3.0
school_event_value 071 0.0% 29.6 29.5 12.5
if_not_fallow 0.070.62 0.0% 0.4 0.4 0.1
value_of_non_garden_land_use 76595 0.0% 353.4 358.2 89.5
costs_of_non_garden_land_use 0.07.3 0.0% 3.0 3.0 1.2
value_of_fallow_non_garden_land 861 0.0% 34.5 34.2 8.8
NPV_garden 9425K 0.0% 2,744.8 2,695.8 777.5
NPV_no_garden 283K 0.0% 751.9 197.1 826.5
decision -1.1K5.0K 0.0% 1,992.9 2,098.4 1,104.6

Summary of the savings

summary(garden_simulation_results$y$decision)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   -1145    1303    2098    1993    2733    5043

Summary of costs

summary(garden_simulation_results$y$total_costs)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   232.2   322.1   343.5   343.8   366.5   447.1

First year

summary(garden_simulation_results$y$Cashflow_garden1)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   65.86  338.74  436.00  444.09  539.97  908.71

Cashflow of the garden option

source("functions/plot_cashflow.R")
plot_cashflow(mcSimulation_object = garden_simulation_results, 
              cashflow_var_name = "Cashflow_garden")

EVPI

# Subset the outputs from the mcSimulation function (y) by selecting the correct variables be sure to run the multi_EVPI only on the variables that the we want.
mcSimulation_table <- data.frame(garden_simulation_results$x, 
                                 garden_simulation_results$y[1:3])
source("functions/multi_EVPI.R")
evpi <- multi_EVPI(mc = mcSimulation_table, first_out_var = "NPV_garden")
## [1] "Processing 3 output variables. This can take some time."
## [1] "Output variable 1 (NPV_garden) completed."
## [1] "Output variable 2 (NPV_no_garden) completed."
## [1] "Output variable 3 (decision) completed."
source("functions/plot_evpi.R")
plot_evpi(evpi, decision_vars = "decision")
## Warning: There are no variables with a positive EVPI. You probably do not need
## a plot for that.

Projection to Latent Structure

Here we run a Projection to Latent Structures (PLS) model, a flexible type of regression model. The results of this model can give us some sense of the correlation strength and direction for model variables and our outcome variables.

source("functions/pls_model.R")
pls_result <- pls_model(object = garden_simulation_results,
                                resultName = names(garden_simulation_results$y)[1], 
                                ncomp = 1)

input_table <- read.csv("inputs_kontum_garden.csv")
source("functions/plot_pls.R")
plot_pls(pls_result, input_table = input_table, threshold = 0)

The full repository can be accessed with the following QR code.